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Abstract 

We propose new summary statistics for intensity-reweighted moment stationary point pro- 
cesses that generalise the well known J-, empty space, and nearest-neighbour distance dis- 
tribution functions, represent them in terms of generating functionals and conditional inten- 
sities, and relate them to the inhomogeneous reduced second moment function. Extensions 
to space time and marked point processes are briefly discussed. 
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1 Introduction 

The analysis of data in the form of a map of (marked) points often starts with the computation 
of summary statistics. Some statistics are based on inter-point distances, others on the 
average number of points in sample regions, or geometric information. For a survey of the 
state of the art and pointers to the literature, the reader is referred to the recent handbook 
of spatial statistics [TO] . 

In the exploratory stage, it is usually assumed that the data constitute a realisation of 
a stationary point process and deviations from a homogeneous Poisson process are studied 
to suggest a suitable model. Although stationarity is a convenient assumption, especially if 
- as is often the case - only a single map is available, in many areas of application, though, 
heterogeneity is present. To account for possible non-stationarity, Baddeley et al. [3] defined 
a reduced second moment function by considering the random measure obtained from the 
mapped point pattern by weighting each observed point according to the (estimated) intensity 
at its location. Gabriel and Diggle [9] took this idea further into the domain of space time 
point processes. 

In this paper, our aim is to define an extension of the J-function |18| that is able to 
accommodate spatial and/or temporal inhomogeneity. The idea underpinning the J-function 
is to compare the point pattern around a typical point in the map to that around an arbitrarily 
chosen origin in space in order to gain insight in the interaction structure of the point process 
that generated the data. The power of the J-function in hypothesis testing was assessed in [6] 
and [26]. Extensions to multivariate point processes were proposed by [19] . and window based 



J- functions suggested by [2] and [6|. For applications in agriculture, astronomy, forestry and 
geology, see [SI HI El [El ED [21] . 

The plan of this paper is as follows. In Section we fix notation and recall some basic 
concepts from stochastic geometry. In Section [3] we describe the most important summary 
statistics that are being used in exploratory analysis of point patterns under the assumption 
of stationarity. Section 0] introduces the new statistic Jmhom and gives representations of it in 
terms of generating functionals and conditional intensities. Section[5]is devoted to the explicit 
computation of Jmhom f° r some important classes of point process models. In Section [6] we 
develop a minus sampling estimator and apply it in Section [7] to simulated examples. The 
paper closes with suggestions for further extensions to space time and marked point processes. 



2 Preliminaries and notation 



Throughout this paper, let X be a simple point process on M. d . Its intensity measure A is 
defined by 



A(B) = E 



.xex 



for Borel sets B C M. d . We assume that A is locally finite, i.e. A(B) < oo whenever B is 
bounded, and absolutely continuous with respect to Lebesgue measure so that 



MB) 



B 



X(x) dx 



for some non-negative measurable function A referred to as intensity function. Heuristically 
speaking A(x) dx is the probability of observing some point in the infinitesimal region dx and 
represents the heterogeneity of X. 

Note that the intensity measure is also known as the first order factorial moment measure 
of X. Higher order factorial moment measures AW, n G N, are defined by 



AW(Bi 



x B n ) = E 



.,x n ex 



l{xi G B x ;... ;x n G B n } 



where the superscript ^ indicates that the sum is taken over all n-tuples of distinct points 
and the Bi are Borel subsets of M d . As the intensity measure, A^ n ) is not necessarily locally 
finite, nor guaranteed to have a Radon-Nikodym derivative. If A^ n ) is absolutely continuous 
with respect to the n-fold product of Lebesgue measures, 



A W (BiX...xB n )= I ••• I 



p (n) (xi, ...,x k )dx 1 --- dx k 



for some non-negative measurable function called n-th order product density of X. Note 
that // n ) is permutation invariant and satisfies the integral equation 



E 



E 



xi,...,x n £X 



g(xi,...,x n ) 



g(xi,. ..,x n ) /? (n) 0i, • • • ,x n ) dx\ 



■ dx r 
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for all non- negative, measurable functions g > 0. Thus, p( n \dxi, . . . , dx n ) may be interpreted 
as the infinitesimal probability of finding points of X at each of dx\, . . . , dx n . For further 
details, see for example [TU 1 H2] , US]. 

In the physics literature, n-point correlation functions tend to be used instead of product 
densities [22J. They are defined recursively by 



6 = i; 



p t - n '(x 1 ,...,x n ) 
X(x l ) ■ ■ ■ X(x n ) 



n 



k=lD u ...,D k 



where the last sum ranges over all partitions {D±, . . . of {1, . . . ,n} in k non-empty, 

disjoint sets, and the X£>. = {xi : i G Dj}, j = 1, . . . , k form the corresponding partition of 
points. Since for a Poisson point process £ n = for n > 1, heuristically speaking n-point 
correlation functions account for the excess due to n-tuples in comparison to a Poisson point 
process with the same intensity function. 



3 Summary statistics 

Summary statistics are used by spatial statisticians as tools for exploratory data analysis, 
testing, and model validation purposes. Popular examples include the nearest neighbour 
distance distribution function G, the empty space function F, the reduced second moment 
function K and the J-function. More specifically, for a stationary point process X with 
intensity A > 0, 

(1) 



F(t) = P(xns(o,t)^0), 

G(t) = P lo (XnB(O,t)^0), 

K(t) = E !0 [£ xex l{z€i?(0,t)}]/A, 

J{t) = (1 - G(t)) / (1 - FQt)) , 



where B(0,t) is the closed ball of radius t > centred at the origin and P !0 denotes the reduced 
Palm distribution of X. For further details about these and other summary statistics, see 
for example |12| . Note that the J-function is defined only for t such that F(t) < 1. Values 
larger than one indicate inhibition, whereas J(t) < 1 suggests clustering, but note the caveats 
against drawing too strong conclusions in [2]. 

All statistics defined in ([1]) can be expressed in terms of product densities. Indeed, if 
the second order factorial moment measure exists as a locally finite measure with Radon- 
Nikodym derivative p^(x\, a^) = /9^ 2 ^(||xi — X2II), 

K(t)= [ ^-Ml d x= I (1 + 6(11x11)) da:. 

JB(0,t) A JB(0,t) 

Clearly, K depends only on product densities up to order two. In contrast, the empty space 
function depends on product densities of all orders [27] . 



^ n\ J Bm J E 



/9 (n) (xi, . . . , X n ) dX! ■ ■ ■ dx r 

B(0,i) JB(0,t) 
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provided all order product densities exist and the series is absolutely convergent, that is, 
hmsup n _ ) . 00 \) ' < 1. Similarly, 



^ 7b(o,*) JB(o,t) A 
provided that the series is absolutely convergent. Thence [17] . 

n=l 

for all i > for which < 1, where J n (t) = f B , o t \ • • • / B ( 0)t ) £n+i(0, xi, . . . , x n ) dx 1 ■ ■ ■ dx n . 
If product densities of all orders do not exist, one may truncate the series. Indeed, using only 
product densities up to second order gives 

J(t)-la-\(K(t)-\B(0,t)\) 

so the iT-function can be seen as a second order approximation to the J-function. 

For non-stationary point processes, the definitions in ([T|) depend on the choice of origin 
and adaptations are called for. To this end, Baddeley et al. [3] introduced the notion of 
second order intensity-reweighted stationarity . A point process X possesses this property if 
the random measure 

is second-order stationary. Here, S x denotes the Dirac measure that places a single point at 
x. Clearly if H is stationary, it is also second-order stationary but the converse does not hold. 
Examples of second order intensity-reweighted stationary point processes include Poisson 
point processes, the random thinning of a stationary point process, and log Gaussian Cox 
processes driven by a Gaussian random field with a translation invariant covariance function. 
Cluster processes, as well as more general superposition processes, typically are not second 
order intensity-reweighted stationary. 

For a second order intensity-reweighted stationary point process, an inhomogeneous K- 
function [3] can be defined by 



^inhom(t) := r^jE 



^ l B (x)l{yeB(x,t)} 

2^x,y&X \(x) X(y) 



regardless of the choice of bounded Borel set B C R rf and using the convention a/0 = 
for a > 0. Indeed, Ki n h ora (t) = tCs(B(0,t) \ {0}), where /C= is the reduced second moment 
measure of the random measure S. 

Gabriel and Diggle [9] restrict themselves to point processes X that are simple and have 
locally finite moment measures of first and second order. Additionally they assume that X 
has an intensity function A that is bounded away from zero and a pair correlation function 

9^y) = 9{\\ x -y\\) = ^^ 



4 



that depends only on \\x — y\\. Clearly in this case 



tinhorn (t) = r4i f f g(\\z\\)dzdx= [ g(\\z\\)dz, 
lBJB(<d,t) JB(0,t) 



\B\ 



which for any inhomogeneous planar Poisson process reduces to nt 2 . 

Baddeley et al. [3j briefly discuss how to define empty space and nearest neighbour 
distance distribution functions for inhomogeneous point processes. First, for given x € R d 
and t > 0, they propose to determine r(x,t) by solving 



B(x,r(x,t)) 



then set 



F x {t) = F{d{x,X) <r(x,t)) 
G x (t) = ¥ lx (d(x,X) <r(x,t)), 

where d(x, X) denotes the shortest distance from x to a point of X. For Poisson processes, the 
above definitions do not depend on x and are both equal to 1 — e~*. The obvious drawback of 
such an approach is that r(x, t) may be hard to compute in practice. Moreover, the definitions 
depend on x as well as t. Our goal in the present paper is to give an alternative definition 
of F, G, and J for intensity-reweighted moment stationary point processes based on their 
representation in terms of product densities that does not depend on the choice of origin and 
is easy to use in practice. 

4 Inhomogeneous J-function 

Let X be a simple point process on M. d whose intensity function A exists and is bounded away 
from zero with inf x A(x) = A > 0. Assume that for all n £ N the n th order factorial moment 
measure exists as a locally finite measure and has a Radon-Nikodym derivative with 
respect to the n-fold product of Lebesgue measure £ with itself for which the corresponding 
n-point correlation function £ n is translation invariant, that is, t; n (xi + a, . . . ,x n + a) = 
£ ra (xi, . . . ,x n ) for almost all a € M. d . We shall call such a point process intensity-reweighted 
moment stationary. Note that a fortiori X is second order intensity-reweighted stationary. 
Moreover, a stationary point process is also intensity-reweighted moment stationary. 

Definition 1. Let X be an intensity-reweighted moment stationary point process. Set 



and define 



J n (t) = •••/ £ n+ x(0,xi,...,x n )dxi---dx n 

JB(0,t) JB(0,t) 

J'mhovaip) = 1 "I - / . j Jnifyi 
n=l 

for all t > for which the series is absolutely convergent, that is, for which limsup ri 
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A few special cases deserve to be mentioned. For a Poisson point process with intensity 
function A(-), as the n-point correlation functions vanish for n > 1, so do the J n whence 
^inhom(^) = 1 f° r all t > 0. Furthermore, if X is stationary, A = A and by [TTJ, Prop. 4.2], 

"tinhorn = J- 

Like in the stationary case considered in Section [31 the series in Definition [T] may be 
truncated, for example when X is only second order intensity-reweighted stationary or not 
all n-point correlation functions exist. For n = 1, we obtain 

Jinhom(i) - 1 » -A / 6(0, a;) da: = -\(K inhom (t) - \B(0,t)\) . 

JB(0,t) 

In the remainder of this section, we rewrite Jinhom m terms of generating functionals and 
conditional intensities. Recall that for any function v : M. d — > [0, 1] that is measurable and 
identically 1 except on some bounded subset of M. d , the generating functional at v is defined 

as 



G(v) = E 



n v(x) 



where by convention an empty product is taken to be 1. The distribution of X is deter- 
mined uniquely by its generating functional [3 Prop. 7. 4. II]. The factorial moment measures, 
provided they exist as locally finite measures, can be derived from the generating functional 
using its Taylor expansion [71 Prop. 7.4. III]. Conversely, if product densities of all orders 
exist, let it be a measurable function with values in [0, 1] that has bounded support. Then 

G(v := 1 - it) = 1 + - j j ■■■ j u(xi) ■ ■ ■ u(x n ) /0 (n) (xi, . . . ,x n ) dx\ ■ ■ ■ dx n , 

n=l U - •> J 

provided the series converges [251 P- 109]. 
Theorem 1. Write, for t > and a € M. d , 

Xl{x G B(a,t)} 



<{x) 



\(x) 



and assume that limsup n _ >00 [K J B(ot) • • • / B ( 0)t ) P \{ x ty'\{xn) dxi ' ' ' dXn ) ^ < lm Under the 
assumptions of Definition^ for almost all a € M. d , 

T. m _G !a (l"<) 

^inhom \") 



for all t > for which the denominator is non-zero, where G' a is the generating functional 
of the reduced Palm distribution P !a at a, G that of P itself. 

Note that for a stationary point process, uf(x) = l{x £ B(a,t)}, hence 
G(l — <) = F(X n B(a, t) = 0) = 1 - F(t). 
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Therefore, the generating functional in the denominator can be interpreted as the inhomo- 
geneous counterpart of the empty space function. A similar interpretation holds for the 
numerator in terms of the nearest neighbour distance distribution function and one retrieves 
the classic definition of the J-function given in Section [3J At this point it should be empha- 
sised that the numerator and denominator in the definiton of Jmhom generalise respectively 
the nearest neighbour distance distribution function and empty space function. 



Proof: We begin by showing that 



A }{Xi € B(x,t)} 

^—<Xi,...,X n £X ll 

1=1 



p( n+1 \0,x 1 ,...x n ) 
(0,t) JB(0,t) H x l) ■ ■ ■ K x n) 



dx\ . . . dx r 



for almost all x G M. d . To see this, consider the functions 

l{x £ A} ^ -A- l{xj £ B(x,t)} 
9A{X > X) = ^ - vll Xjx-) 

defined for all bounded Borel sets A C M. d . By the definition of Palm distributions and the 
Campbell-Mecke formula, 



g{x,tp) \{x) dP lx (p) dx = E 



^Tg(x,X\{x} 
.xex 



Using Fubini's theorem, for our choice of g the left hand side can be written as 

rr l{xi € B(x,t)} 



E 



E 



xi,...,x n £X 



n 

i=\ 



X(xi) 



da- 



while the right hand side is equal to 

E 



^# l{x € A} yj l{ X j € B(x,t)} 

\(<r\ 11 



X(x) 1J - X(x 
i=i 

The expectation can be computed in terms of p^ n+l ^ and equals 

p (n+1) (x,xi, ...x n 



I \ JB{x,t) JB(x,t) X(x) X{xi) ■ ■ ■ X(x n ) 

p( n+1 \0,x u ...x n ) 




dx dx\ ■ ■ ■ dx n = 
dx dxi ■ ■ ■ dx n 



i A JB(o,t) JB(o,t) X(0) X{xi) ■ ■ ■ X(x n ) 
by the translation invariance of the n-point correlation functions. Hence 



E 



xi,...,x n £X 



n 

i=i 



i{ Xi e B(x,t)} 
X(xi) 
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is constant for almost all x E 
Next, note that 



n / \l{x £ B(a,t)} \ ~ (_A)*W " l{ Xl £B(a,t)} 

H V A(x) / ^ n! ^ a! i I ...,xn6xll 



A(xi 



Since the number of points in X D -B(a, i) is almost surely finite, the expressions are well- 
defined under the convention that an empty product takes the value one. Consequently, for 
almost all a, 



(2) ^(i-o = i+f;^/ •••/ 



p^+ 1 \0,x 1 ,...,x n ) 
{ nl is(o,t) /fl(0,t)^)%)-%) 

provided the power series in the right hand side is absolutely convergent. 
By the discussion preceeding the statement of the theorem, 



dx\ ■ ■ ■ dx r 



(3) 



^ n! y B (o,t) ./£ 



/? (n) (xi, . . . ,x ra ) 
B(o,t) 7s(o,t) A(xi) • • • A(x n ) 



dxi ■ ■ ■ dx n , 



since the power series in the right hand side is assumed to be absolutely convergent. 

Upon recalling the definition of the n-point correlation functions and splitting into terms 
that do or do not contain the origin, one obtains that the right hand side of (|2|) is equal to 



i+E^r E •«*> E E 

n=l DC{i,...,n} fc=l Dl ... Dfc#0 disjoint 

UDj={l,...,n}\£) 



4(D X ) • • • 4(£>fc) 



(with E^ =1 = 1) which in turn can be written as 



(4) 



where 



i+E 



n=l 



nl 



i+E 



E E 

ml ^ z — ' 

fc=1 o 1 ,...,r> fe ^0 disjoint 

UDj={l,...,m} 



^ n (xi, . . . ,x n ) ctei • • • dx r . 



4(-Di) • ■■ I n{D k ) 



>B(0,t) JB(0,t) 

and n(-D) denotes the cardinality of the set D. The sum over k in the rightmost term of 
can be written as 



p^(x l ,...,x m ) 



'B{a,t) JB(p,t) A(xi) ■ ■ ■ A(x m ) 
hence the second term in ((4]) is equal to the right hand side of (|3]). Finally, since both sums 



in (Hj) are absolutely convergent, so is (|2|), an observation that completes the proof. 



□ 
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Next, we focus our attention on conditional intensities X(x;X), x € M. d . Assuming they 
exist, they are defined in integral terms by 



E 



5>(*,X\{*}) 



xex 



K- x [g(x, X)} X(x) dx = / E [g(x, X)X(x; X)} dx 



for any non-negative measurable function g. 

Theorem 2. Assume that X admits a conditional intensity and define the random variable 
W a (X) := IXrex (1 ~ u t( x ))- Then, under the assumptions of Theorem^ E [W a (X)] = 
implies E [X(a;X) W a (X) / X(a)] = for almost all a € M. d , and otherwise 



^inhom \P) — E 



A (a; X) 



W a (X) 



/EW a (X), 



the W a -weighted expectation of A(a; X) j A(a) . 



Consequently, Jmhom(i) < 1 Cov ^ A ^ a ^ , W a (X) ) < with a similar statement for the 
opposite inequality sign. 

Proof: Consider the functions 



9a(x,X) 



1{x £ A} 



n 1 



yex 



Xl{yeB(x,t)} 



defined for all bounded Borel sets A C M. d . Arguing as in the proof of Theorem Q] and using 
the definition of conditional intensities, one obtains 



f E te TT A Alfa £ B(x,t)} \ 
J A y \ l x \ A(y) J 



dx 



E 



X(x; X ) 



n(i 



Hence, 



E ! ' 



yex 



Xl{yeB(x,t)} 



E 



A(x;X) j-j- /_ Xl{yeB(x,t)} 



X(x) 



yex 



A(y) 



for almost all x G An appeal to Theorem [T] completes the proof. 



□ 



5 Theoretical examples 
5.1 Poisson process 

Let X be a Poisson point process with intensity function A : M. d — > R + that is bounded away 
from zero. Since p^ n \x±, . . . ,x n ) = rXjA(xj), the n-point correlation functions vanish for 
n > 1, SO Jinhom(^) = 1 for all t > 0. 



9 



The generating functional of X is 



G(l — u) = exp 



u(x) X(x) dx 



In particular, for the function u = defined in Theorem [Tj G(l - u° t ) = exp [-A|B(0, t)\] . 
Also, since according to Slivnyak's theorem for a Poisson point process P !0 = P, G'°(l — u^) = 
G(l — u®). Finally, the conditional intensity A(-,X) of X coincides almost everywhere with 
the intensity function A(-). 

5.2 Location dependent thinning 

Let X be a simple, stationary point process on W 1 for which product densities p^ of all 
orders exist. Let p : M d — > (0, 1) be a measurable function that is bounded away from zero 
and consider the thinning of X with retention probability p(x) as in Example 8.2 of [7\. Since 

(n) 

the process is simple, the product densities p th of the thinned point process can be expressed 

in terms of those of X by p{^{xi, . . . , x n ) = p^ n \x\, . . . , x n ) Yli=iP(. x i)- ^ n particular, the 
intensity function of the thinned point process is Ath(a^) = Ap(x), where A > is the intensity 
of X. Consequently, 

p^ } (xi,...,x n ) = p {n) {x 1 , ...,x n ) 

AthOl) • • • AthOn) A n 

Therefore, the n-point correlation functions of the thinned point process coincide with those 
of the underlying stationary point process X, £* h (xi, . . . , x n ) = Cni^i, ■ ■ ■ , x n ), and inherit the 
property of translation invariance. Hence J^{t) is equal to the J n -function of the underlying 
point process X. As the intensity function of the thinned point process is bounded from 
below by \p where p is the infimum of the retention probabilities, 



oo 



^bihomW — 1 + ^ " — — p Jn{t) 

n. 

n=l 

for all t > for which the series converges. Note that the power series coefficients are identical 
to those in the power series expansion of the J- function of X. 

The generating functional of the thinned point process is Gth{v) = G(vp + 1 — p), where 
G is the generating functional of X. Hence 

G th (l - 4) e B(0, t)} \ =G(l-p 1{- e B(0, t)}) = E [(1 - p)"(*nB(o,t)) 

the generating function of the number of points of X that fall in B(0, t) evaluated at 1 — p. 

As the reduced Palm distribution of the thinned point process coincides with a random 
location dependent thinning of the reduced Palm distribution of X with retention probabilities 
given by the function p, 



<(l-n t °)=E !0 



(1-p) 



n(xns(o,t)) 
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so that under the assumptions of Theorem [TJ 



K 10 [(1 _ p) n(XnB m) 

e Anhom\'v 



E [(1 -p)n(XnB(0 ; 



0)1 



To conclude this example, note that the assumption of stationarity of the underlying point 
process X may be weakened to intensity-reweighted moment stationarity. 



5.3 Scaling 

Let X be a simple point process on W 1 for which product densities p( n ' of all orders exist. 
Let c > be a scalar constant and map the point pattern X to cX. Then all order product 
densities p£g of cX exist and are given by p^(x%, . . . , x n ) = c~ dn p( n \xi/c, . . . ,x n /c). In 
particular for n = 1, X c x{x) = c~ d X(x/c). Therefore the n-point correlation functions 

£n 0^1 j ■ ■ ■ j "Efi) = £n (2<l/c, . . . , X n j c) 

of cX are invariant under translations if and only if the n-point correlation functions £ n of X 
are, in which case the J„-functions J^ x of cX are scaled versions J^ x (t) = c dn J n (t/c) of the 
corresponding functions of X. Furthermore, inf x€R d X c x(x) = Xc~ d , so the inhomogeneous 
J-function of cX is 

°° (— \r~ d \ n °° (— \) n 

JELnW = 1 + E , Cdnj "( t / C ) = 1 + E ^^T" W<0 = Jinhc m (t/C), 

z — ' n! z — ' n! 

n=l n=l 

the inhomogeneous J-function of X evaluated at t/c provided the series converges. Note that 
in contrast to the thinning case, the power series coefficients are not identical to those of the 
underlying point process X. 

The generating functional of the scaled process is given by G c x(v) = G(v(c-)), where G 
is the generating functional of X, whence 

Gcx (l - c J A( A /c) l{- € B(0,t)}) =c(l- e B{Q,t/c)}) . 

Similarly, noting that d¥'f x ((p) = dF- x / c ((p/c), 

G'tx (l " c _^ A /c) l{- e fl(a,t)}) = G !a / C (l - iylf € B(o/c,t/c)}) . 

To conclude this example, a conditional intensity of is obtained by scaling that of X, i.e. 
A c x(a;,v 9 ) = c~ d X(x/c,ip/c) [11], from which we retrieve the formula </^h om (i) = ^mhom (*/ c ) 
under the assumptions of Theorem [2j 
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5.4 Log Gaussian Cox process 



Write Q for the distribution of a random measure defined in terms of its Radon-Nikodym 
derivative A with respect to Lebesgue measure. We assume that all moment measures of 
the random measure exist and are locally finite. Let X be the Cox process directed by the 
random intensity process A, that is, given a realisation A = A, X is a Poisson point process 
with intensity function A. It follows from [TJ p. 262] that the factorial moment measures of 
X exist and are equal to the moment measures of the driving random measure. Hence X has 
product densities p( n \xi, . . . , x n ) = E JlILi &-{ x i)] ■ Moreover, the reduced Palm distribution 
of X at x is the distribution of a Cox process with driving random measure distributed as 
Q x , the Palm distribution of the driving measure of X at x \25\ p. 141]. 

The class of log-Gaussian Cox processes [20J is especially convenient. For models in this 
class, 

A(x) = exp [Z(x)] 

where Z is a Gaussian field. Such a field is defined fully by its mean and covariance function. 
Write p(x) for the mean function, cr 2 (x) for the variance of Z(x) and r(x,y) for the corre- 
lation function. In other words, the covariance function of Z is given by a(x) a(y) r(x,y). 
Conditions have to be imposed on these functions in order to make the resulting Cox pro- 
cess well-defined. In particular, the intensity function must be integrable almost surely, and 
^a(B) = Jg A(x) dx a finite random variable for all bounded Borel sets B C W 1 . Moreover, 
the distribution of the random measure ^\ must be uniquely determined by that of Z. Suf- 
ficient conditions are given in [TJ Thm. 3.4.1] for zero mean Gaussian processes. Therefore, 
we additionally assume that the mean function \i is continuous and bounded. Now, since 
E A-(xi)] = E [e^ j , the moment generating function of the normally distributed 
random variable ^2 i Z(xi) evaluated at 1, 



/9 (n) (a;i, ...,x n ) = exp 



Yl [p( x i) + 2 J +^2^(x i )a(x j )r(x i ,x j ) 

i=l ^ ' i<j 



Specialising to n = 1, it follows that the intensity function is logA(x) = + a 2 (x)/2 
whence 

p( n \xi, ...,x n 



X(xi) ■ ■ ■ \{x n ) 



exp 



^cr(xj) a(xj) r(xi, : 



i<j 



Thus, if it(-) = g > and r(x, y) = r(x — y), X is intensity-reweighted moment stationary and 
the intensity function is bounded away from zero with infimum exp [a 2 /2 + inf xgR d n(x)] . 

In order to derive an explicit formula for Jmhonn we turn to the generating functional. 
Recall that a Cox process has a generating functional Prop. 8.5.1] defined by G(v) = 
Eg exp [— J(l — v(x)) A(x) dx] . Therefore, for the log-Gaussian Cox process 



G(l-u° t ) =E z exp 
where ft denotes mf xGR d e^ x >. 



-H / e z( - x) ~^dx 

'B(0,t) 
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The Palm distributions Q x of a log Gaussian random measure are A(x) = e z ( x )-weighted. 
To see this, note that the Campbell measure evaluated at a bounded Borel set B (ZM. d and 
F in the smallest cr-algebra for which ^\(B) = j B A(x)dx is finite for all such B, can be 
calculated as 



C(B xF)=E Q [1 f (* a) tf A (B)] = I \{x) \ [ 1 F (* A ) ^ d Q(* 

Jb u K x ) 



dx 



by Fubini and the existence of a u-finite intensity measure A(-) that is bounded away from 
zero. Therefore, 



G !a (l-<) 



■ exp 



-\x I e 

I B(a,t) 



dy 



Since Y(x) := Z(x) — fJ<(x), x S R d , is a stationary Gaussian process, the above generat- 
ing functional does not depend on the choice of a. Therefore, under the assumptions of 
Theorem [TJ 



J; 



inhom 



(t) 





exp 




^/B(o,t) e Y(x) dx 






Ey [e y (°] Ey exp 


-A*/B(o,t) e Y(x) dx 



The mixed Poisson process considered in [18J is a special case. 

Note that Jmhom(^) < 1 if and only if the random variables e Y ^ and e H ' JB (o,t)' are 
negatively correlated. The geostatistical models used in practice, for example the one we 
shall use in Section [TJ have a positive, continuously decreasing, correlation function. There- 
fore, by Pitt's theorem [23], the Gaussian fields defined by such correlation functions are 
associated. Under the further conditions of [1, Thm. 3.4.1.], the sample functions Y(-) 
are almost surely continuous and hence the integral of e Y over B(0,t) is uniquely defined 
and the limit of Riemann sums over ever finer partitions of B(0,t). Since Y is associated, 

Gov (e y (°\ e" c » £< eV( " l) ) < for all finite sums with positive scalar multipliers Cj > 0. Upon 
taking the limit, it follows that Ji n hom(^) < 1- 



6 Estimation 

The goal of this section is to develop an estimator for the inhomogeneous J-function of Defi- 
nition [TJ For this purpose, we shall use the representation in terms of generating functionals 
of Theorem [TJ and apply the minus sampling principle outlined in [251 p. 127]. 

Specifically, let Let W C M d be a compact set with non-empty interior and suppose the 
point process X is observed in W. For clarity of exposition, we assume that the intensity 
function A is known. If it is not, it can be estimated (for instance using kernel estimation [5]) 
and plugged into the estimators outlined below. 

In order to estimate G(l — u®), let L C W be a finite point grid. Set 



(5) G{\-u 



J2i k eLnw et T\xexnB(i k ,t) 
1 ' := #L n W et 



A 

A(as) 
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where W et is the eroded set {x G W : d(x,dW) > t} = {x G W : x + B{0,t) C W}. 

Note that G(l — u®) is an estimator as for all grid points Ik £ Wet the ball B(lk,t) is fully 
contained in W so that no points of X \ W are needed for the computation of the product in 
the numerator of ([5]). 
Similarly, set 



(6) 



0(l-<) 



x k exnw e t T[xex\{x k }nB(x k ,t) 



1 



A 

A^y 



#x n w et 

Compared to ©, the grid points l/. are replaced by the points x^ of X n Wet- Again, 
a function of X n W only. 



is 



Proposition 1. Under the assumptions of TheoremUl the estimator |5]j is unbiased, is 
ratio-unbiased. 



Proof: We claim that 
(7) E 



n 

xexnB(i k ,t) 



1 



X(x) 



G(l - «?) 



for all 4 £Ln Wet. To see this, note that 

A 



«» n 

xGXns(z fc ,t) 



1 



A(x) 



n=l 



n 

i=l 



l{zj-jfc £ £(0,t)} 
A(xi) 



Hence, under the assumptions made, the expectation of ([7]) can be expressed as 



(-A)" 



n=l 



tt! 



B(i fc) t) 



/? (n) (^l, 



B(Z fc ,t) 



nr=i A(xi) 



dxi • • • dx r 



which, because of the translation invariance of the integrands, reduces to 



(-A)' 



n=l 



tt! 



B(0,t) 



p( n )(q:i,...,a: n ) 
B(o,i) n™=i Kxi) 



dx\ ■ ■ ■ dx n = G(l — u®). 



This proves the claim, from which the unbiasedness of © follows. 

Next, turn to the numerator of ©. By the definition of Palm distributions and the 
reduced Campbell-Mecke theorem [25j p. 107], its expectation can be expressed as 



/ / a (-) n f 1 

J Jwet yetp \ 



By ([8]), Fubini, the first equation in the proof of Theorem [TJ and (|2j), the Palm expectation 
in the integrand equals G !0 (l — ) for almost all x, hence the numerator of © equals 
G !0 (l — u®) J WQt X(x) dx. As the expectation of the denominator is equal to J WQt X(x) dx, 
(El) is ratio-unbiased as claimed. □ 
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Figure 1: Each row contains a realisations of a point process in the leftmost frame, the graphs 
of © (solid line) and ([6]) (dashed line) in the middle frame, and the graph of ^mhom(^) (solid 
line) compared to irt 2 (dashed line) in the rightmost frame. The models are a Poisson point 
process (top row) , a log Gaussian Cox process (middle row) , and a thinned hard core process 
(bottom row). 
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7 Examples 



In order to see how Jmhom(i) works in practice, we simulated realisations of three of the 
models presented in Section [5j Typical patterns are displayed in the leftmost column of 
Figured! In all three images a smooth intensity gradient can be observed: more points are 
located near the bottom of the square than near the top. However, the interaction structure 
seems different. For example, the middle picture contains groups of points that are close 
together, with large gaps in between the clusters. In the lower picture on the other hand, 
points seem to avoid being very close together and are more evenly spaced out. In the top 
picture, both very small and very large interpoint distances occur. In order to quantify the 
above qualitative remarks, we applied the ideas presented in this paper and compared the 
results to those obtained by a second order analysis. To simulate the patterns and calculate 
the estimators, the R packages spatstat0 and RandomField^ were used. 

Poisson point process The first example is a heterogeneous Poisson point process with 
intensity function A (x, y) = 100 e~ y . Note that the mean number of points is 100(1 — e _1 ) ~ 60 
per unit area. A realisation is shown in the top left frame in Figure [TJ The top middle frame 
shows © (solid line) and ([6]) (dashed line). It can be seen that the graphs lie close together, 
in accordance with the fact that for any Poisson point process, Jmhom = 1- For comparison, 
the plug in minus sampling estimator of i^mhom is shown as the solid line in the top right 
frame. Again, the graph is close to that of the theoretical value 7r t 2 (dashed line in the top 
right frame). 

Log Gaussian Cox process The second example is a log Gaussian Cox process. The 
defining Gaussian random field has exponentially decaying correlation function, unit variance, 
and mean function ji satisfying e^ x,y ' = 100 e -y-1 / 2 . Note that the intensity function of the 
Cox process thus defined coincides with that of the Poisson point process discussed above. 
A realisation is shown in the middle row's leftmost frame in Figure [TJ The middle frame 
in the same row show ([5]) (solid line) and ([6]) (dashed line). Note that the graph of (|6|) lies 
well below that of ©, indicative of attraction between points due to the positive correlation 
of Z after accounting for the inhomogeneity. For comparison, the plug in minus sampling 
estimator of i^inhom is shown as the solid line in the rightmost frame in the middle row. From 
about t = 0.13, the estimated value is smaller than irt 2 . 

Thinned hard core process The third example is a thinned hard core (Strauss) process 
defined by its conditional intensity /3 l{d(x, X\{x} > R}. A realisation for (3 = 200, R = 0.05 
and retention probability p(x, y) = e~ y is shown in the bottom left frame in Figure [TJ The 
middle frame in the bottom row show ([5]) (solid line) and (|6|) (dashed line). Note that the 
hard core distance is clearly reflected in the flat initial segment in the graph of ([6]), which lies 
abovethe graph of ([5]) up to about r = 0.2, indicative of the inhibition between points due to 
that present in the underlying hard core process after accounting for the inhomogeneity. For 

1 Adrian Baddeley, email: adrian@maths.uwa.edu.au; and Rolf Turner, email: r.turner@auckland.ac.nz 
2 Martin Schlather, email: martin. schlather® math.uni-goettingen.de 
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comparison, the plug in minus sampling estimator of i^mhom is shown as the solid line in the 
bottom right frame. The estimated value is smaller than that of a Poisson point process up 
to about t = 0.2 confirming the picture painted by the Ji n h om -function approach. 



8 Summary and extensions 

In this paper, we defined a J-function for intensity-reweighted moment stationary point 
processes in terms of their n-point correlation functions and gave representations in terms 
of the generating functional and conditional intensity. We calculated Jmhom explicitly for 
the three representative classes of intensity-reweighted moment stationary point processes 
presented in |3j, derived an estimator, and presented simulation examples. 

Although this paper focussed on point processes on R rf , the approach may be extended 
to space time or marked point processes. First, assume that Y is a simple point process on 
the product space M. d x M equipped with the supremum distance whose intensity function 
A(-) exists and inf^) \(x,t) > 0. Furthermore assume all order factorial moment measures 
exist as locally finite measures that have Radon-Nikodym derivatives with respect to 
the n-fold product measure of £ with itself, n £ N, and the corresponding n-point correlation 
functions are translation invariant in both components. Define J n as in Definition [TJ from 
which an inhomogeneous space time version of the J-function can be defined. If the series is 
truncated at n = 1, one obtains 



which corresponds to the approach of Gabriel and Diggle [9]. If space and time are 
scaled differently, see Section [573^ J-^ om (t,s) becomes a function of two variables, one for 
spatial distances, the other for time differences, which is more natural in many applications. 

For marked point processes on M. d with marks in some Polish space M equipped with a 
finite reference measure make the same assumptions as above for space time point processes 
except that the n-point correlation functions are required to be translation invariant in the 
spatial component only. For any Borel set B C M and n G N, set Jjf (t) equal to the common 
value of 






which can be estimated usin; 



■g 



minus sampling ideas. 
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